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"^ . Abstract 
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1 Introduction 



H ' Leon Ehrenpreis was a mathematician of remarkable strength, famous 

accomplishments, and extremely wide interests. In the last couple of 
decades of his life, integral geometry was one of the main areas he 
was interested in. This had lead in particular to his involvement with 
problems of tomography [71 [53]. E. T. Quinto and the author had 
the honor of writing at Leon's request an Appendix [82] dedicated 
to tomography for his last book [52]. He was also interested in new 
developments in medical imaging, which in some instances turned out 
to be directly related to some integral geometric and PDE problems 
he considered in [52] . It is thus appropriate to address some of these 
new techniques in an article dedicated to Leon's memory. 
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It is natural to wonder, why do we need new methods of medical 
imaging in the first place, if we already have the whole bunch of well 
developed ones |94 1ll03j ? Indeed, one can mention the widely known 
"standard" X-ray CT (Computed Tomography) scanners, MRI, and 
ultrasound scanners, which can be found in most hospitals nowadays. 
There are also less known to the general public, but well developed by 
now PET (Positron Emission Tomography) and SPECT (Single Pho- 
ton Emission Computed Tomography) techniques. Optical Tomogra- 
phy (OT), Electrical Impedance Tomography (EIT), Elastography, as 
well as quite a few others. 

It seems that even the existing methods are way too many. Why 
would one need all of them? One of the main answers is that different 
imaging techniques "see" different things. One can say (in a crude 
approximation) that the CT scan can distinguish between the tissues 
of different density. What if two tissues (e.g., the cancerous and the 
healthy one) have essentially the same density, but absorb significantly 
different amounts of light in certain part of the spectrum? One can 
think of using OT then rather than the X-ray CT. Some modalities 
can show the metabolism (e.g., the PET and SPECT), while some 
others cannot. Some new methods can show the level of oxygenation 
in blood, while those relying upon density would not be able to do 
so. This list can go on and on, and so the reader can see the point in 
having a variety of imaging techniques. 

There are several other parameters that make a difference when 
using different types of scanners. Here are the most common ones: 

1. Safety for the patient and practitioner. Indeed, X-ray scans 
are clearly not too safe, while for instance OT or ultrasound 
tomography are not harmful. If there is a choice, one surely 
would shoot for a safer method. 

2. Cost, in the times of the high medical expenses, is clearly one 
of the major issues. Some imaging techniques (X-ray CT-scan, 
MRI, PET, and some others) require very expensive devices, 
with the price tags in millions of dollars. Some others, however, 
e.g. OT and EIT, are orders of magnitude cheaper. 

3. Contrast is another important feature. For instance, if the pa- 
rameter that a method can detect is, say, electric conductivity of 
the tissue, then one will be able to distinguish between the tissues 
that have a significant conductivity contrast and will not see any 
difference between tissues that happen to have very close electric 



conductivities. So, high contrast between the features that we 
would hke to distinguish, is crucial. One clearly would prefer the 
contrast that is on the order of hundreds (or at least dozens) of 
percents, while one percent contrast, albeit still usable, is much 
less desirable. For instance, some breast tumors on early stages 
might have almost no contrast with the healthy tissues with re- 
spect to ultrasound propagation, but a huge contrast (several 
hundreds of percents) in their optical and electric properties. 

4. Resolution, which determines what is the smallest size of a fea- 
ture that a method can "see," is another very important param- 
eter. Indeed, resolution of several centimeters probably is not 
good for early breast tumor diagnostics, where sub-millimeter 
resolution is desirable. 

These are the reasons, why the quest for new and "better" (at 
least in some of the parameters) imaging methods not only does not 
subside, but intensifies in the recent decades, involving new physics, 
engineering, and mathematics ideas. However, all attempts to find a 
"magic wand" method that would "see everything" and feature low 
cost, safety, high contrast and high resolution, are clearly futile. Thus, 
having a variety of techniques at our disposal is apparently the way 
to go. 

For a mathematician, however, the main motivation for working 
on these various tomographic modalities is that they bring about a 
large variety of challenging and beautiful mathematical problems that 
involve more or less all areas of mathematics. 

Let us describe now the structure of this article. After a gen- 
eral discussion of hybrid methods in Section [2l we will concentrate in 
more details on three of the emerging hybrid methods. The largest 
Section [3] contains an overview of the most developed (both experi- 
mentally and mathematically) Thermo-/Photo-acoustic Tomography 
(TAT/PAT). The next, shorter, Section H] addresses the much newer 
technique of Acousto-Electric Tomograph}{^. Although AET was ini- 
tially suggested (and its principle experimentally proven) by biomed- 
ical engineers [149j . here mathematics is getting developed fast in the 
recent few years, while experimentalist still struggle with reaching 
good signal-to-noise ratios (SNR) in the measurements. This situ- 
ation is reversed in the Ultrasound Modulated Optical Tomography 



^Besides the name AET, suggested in the original paper |149| . other names are also 
used: Ultrasound Modulated EIT (UMEIT) [17], Impcdiography [10], and some others. 



(UMOT, also called UOT), which by now is significantly studied ex- 
perimentally, while the mathematics of this modality is still in its 
infancy (and even mathematical models are still being agreed upon). 
This topic is addressed in ever shorter Section [6l 

Both AET and UMOT rely upon an assumption of a "perfect focus- 
ing" of ultrasound at a given location, which is a crude approximation 
to the reality (see, e.g. [66] )• The synthetic focusing technique, dis- 
cussed in Section [SI allows one to use more realistic sets of ultrasound 
waves. 

The common feature of AET and UMOT (as well as some other hy- 
brid imaging techniques, such as GDI and MREIT, which are just men- 
tioned in this text) is that the measurements provide the researcher 
with some interior information, i.e. a function of an interior location 
X. There has been a rather common feeling that such an interior in- 
formation might stabilize the notoriously unstable modalities such as 
EIT or OT. This issue is briefly discussed in Section [T[ 

The topics surveyed in this article are highly technical and involve 
a wide range of interesting mathematics. Due to space limitations, the 
author was forced to show very few technical details and instead to try 
to give a hand-waving heuristic description. The literature references 
(as well as the references therein) will provide the interested reader 
with more details. 

2 Hybrid imaging methods 

As we have indicated in the introduction, each of the available imaging 
methods has its advantages and deficiencies. For instance, in breast 
imaging ultrasound provides a high (sub-millimeter) resolution, while 
suffers from a low contrast. On the other hand, many tumors absorb 
much more energy of electromagnetic waves (in some specific energy 
bands) than healthy cells. Thus, using such electromagnetic waves 
offers very high contrast. Alas, the resolution in this case is very low. 
One can go on and on with such examples. 

Since both the advantages and disadvantages of various modalities 
come from the underlying physics, can one do anything to combine the 
advantages and simultaneously alleviate the problems associated with 
different type of physical waves/radiation involved? The natural idea 
is to try to this end to combine different imaging modalities into some 
kind of "hybrid" ones. This is how the hybrid (also called multi- 



physics or multi-wave) techniques have been appearing in tiie last 
decade - decade and a half [T0l[T6l[T9l [T2n[T32l[T35l[T39lll4Ulll42j . 

The "hybridization" can occur at different stages of the imaging 
process. Let us recall the crude scheme of all CT methods: on the 
first step, some wave(s) are sent through the body and the outgoing 
(transmitted or reflected) waves are measured; on the second stage, 
mathematical processing of the measured data is done; finally, the 
third stage provides a picture (tomogram). Correspondingly, one can 
combine different techniques at different stages. The simplest, and 
already industrially implemented (e.g., in PET/CT scanners) idea is 
to run two types of scans of the same patient and then somehow 
"combine" the images so they hopefully complement each other. Here 
reconstruction of both individual images requires neither new types of 
scanners, nor new mathematical reconstruction algorithms. Certainly, 
some (often non-trivial) processing is needed to correctly overlap the 
two images (the so called image registration). 

Another option is to combine the techniques on the second stage. 
Namely, after collecting data from two independent types of scans, 
a reconstruction algorithm is applied that uses both data sets. The 
additional information can significantly improve the quality (stability, 
resolution, etc.) of the resulting picture (tomogram). Such combined 
techniques might not require any new physics or engineering, but 
do demand new mathematical processing. There are several recent 
imaging procedures that successfully implement this idea. Some of 
them (CDI - current density imaging pHTIIn] and MREIT [T22l[T43] ) 
combine MRI and EIT measurements. Having the extra MRI data 
makes the mathematical problem of EIT reconstruction (the so called 
Calderon problem, or inverse conductivity problem |15 1 [27 1[3H I37 1 I4H 
[73l[Ml[971 [T30l[T33lll34lll36] ). known for its severe ill-posedness, signif- 
icantly less ill-posed and thus allows for good quality reconstructions 
of the internal electrical conductivity maps. Another actively devel- 
oping method of this kind, called MRE |95ll96j . combines MRI with 
elastography: MRI allows one to observe propagation of elastic waves 
through the tissue, which then leads to mathematical reconstruction 
of mechanical properties (e.g., stiffness) that carry important medi- 
cal diagnostic information. Probably the oldest such combination is 
of CT and SPECT. The CT scan provides the reconstruction of the 
attenuation map, which is used then to recover the distribution of a 
radio-pharmaceutical inside the patient's body |33j . 

Due to the author's limited expertise and lack of space, these types 



of hybrid methods will not be discussed in the article. The interested 
readers are referred to the literature cited above. 

We reserve in this text the name "hybrid methods" only for the 
techniques that combine different types of waves already on the first, 
scanning stage. In the examples that we will discuss this will lead 
to one physical type of irradiation triggering or modulating the other 
one and thus producing new types of measurements, which hopefully 
allow one to irnprove images in comparison with the two techniques 
done separateljo. 

We now move to considerations of some of the hybrid techniques 
in more detail. We start in the next section with the probably best de- 
veloped (both experimentally and mathematically) among the hybrid 
techniques Thermo-/Photo-acoustic tomography (TAT/PAT) and then 
move to the less studied ultrasound modulated electrical and optical 
tomography. 

3 TherrQO-/Photo-acoustic tomography 
(TAT/PAT) 

As we have already mentioned, in many medical diagnostics situ- 
ations, ultrasound displays low contrast (and thus sees the tissues 
as almost homogeneous), while providing fine resolution. Optical or 
radio-frequency EM illumination, on the other hand, gives an enor- 
mous contrast between the cancerous and healthy cells, while both 
are known to suffer from low resolution. How can one combine their 
strengths? The answer is in the photo-acoustic effect, which was dis- 
covered by Alexander Graham Bell [28], but had to wait for another 
century for its applications to follow [32] . 

Imagine that a biological object is irradiated by a wide, homoge- 
neous, but extremely short electro-magnetic pulse in radio-frequency 
range (Fig. [T]). Some part of the electromagnetic (EM) energy will be 
absorbed throughout the tissues. Let us denote by f{x) the density 
of energy absorbed at a location x. It is known that the values of this 
absorption function will be several times higher at the cancerous loca- 



tions than in the surrounding healthy tissues [5nl[88ll8mfin8] lin9 t ll31 



^The reader should be aware that this classification of hybrid modalities into three 
classes, although being reasonable, is not commonly accepted and is used in this text only 
for the author's convenience. 
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Figure 1: The TAT procedure. 



ll32yi39Ml42pi44] . Thus, if this function were known, it would provide 
an extremely valuable diagnostic information. However, the radio fre- 
quency waves are too long to lead to any reasonable resolution. This 
is where the photo-acoustic effect kicks in: the EM energy absorption 
leads to heating, and higher levels of absorption lead to more heating. 
In turn, the resulting thermoelastic expansion creates a pressure wave 
p{x,t) (acoustic wave), whose initial value is essentially (proportional 
to) the function of interest f{x). Placing an ultrasonic transducer (a 
microphone) at a location y at the boundary of the object, one can 
measure the value of p(y, t) at this point for any time t > 0. If we 
now surround the object by transducers located on an observation 
surface S (this can be done using optical interferometers, which do 
not obstruct the acoustic wave propagation), we can collect the values 
oip{y,t) for all {y,t) € 5 x M+ (Fig. [1]). Now the problem reduces to 
finding the initial values p{x, 0) inside the volume surrounded by the 
surface S. This is the idea of Thermoacoustic Tomography (TAT), 
which found its implementation in the middle of 1990s |75 1ll081ll09j . 
in particular in the work of R. Kruger [75], who started a company 
OptoSonics manufacturing TAT devices. The Photo-acoustic version 
(PAT) differs by the choice of heating radiation - a laser beam in- 
stead of radio- waves. A large part of the corresponding mathematics 
is parallel in TAT and PAT cases, and thus we will only mention TAT 
here. 

The mathematics of TAT/PAT reconstruction happens to be fas- 
cinating and occupied attention of a large group of mathematicians 



throughout the 1st decade of this century. The reader can consult 
with the recent surveys, coUections, and books [3l[T2 t [T9 l [56 l [58 l [77 y80^ 
[83l[IIIl[IIll[l2ll[l2ll[l32l[l35l[l39] for details and further references. 

We will now describe the mathematical model of TAT. Let us de- 
note by c{x) the sound speed in the interior of the body. Then the 
pressure p satisfies the following wave equation: 

'gp=c2(x)Ap, inM3xM+ 
piO,x) = f{x), (1) 

lt(0,x)=0 

The data g measured by a TAT machine provides the values of the 
pressure p on the observation surface S: 

f{x) ^ g{y,t) ■.= p \sxR+ ■ (2) 

The goal of TAT thus is inverting this forward operator. 

3.1 TAT/PAT and (restricted) spherical mean 
operators 

In the case of an acoustically homogeneous medium (i.e., when c{x) = const), 
one can reduce the TAT problem to an equivalent integral geome- 
try question. Namely, using the standard Kirchhoff-Poisson formu- 
las |46U71| for the solution of ([T|), the TAT inversion can be reduced to 
the equivalent (see [6]) problem of recovery the function f{x) from its 
spherical averages over spheres centered at the transducers' locations, 
i.e., on the observation surface S: 

f ^ g{y, t) := {Rsf){y, t) = f f{x)da{x),y € 5, t > 0. (3) 

J\x—y\=t 

Here Rs denotes the operator of taking spherical averages over all 
spheres centered on S. 

One can notice that such transforms were considered in more gen- 
eral situation and without any relation to tomography by Leon Ehren- 
preis in his book |52j and by V. Lin and A. Pinkus in |9ip92j for the 
needs of approximation theory and neural networks. The restricted 
spherical mean transforms also play important role in the Radar and 
Sonar studies |40y93|. One can find a list of other applications in [6]. 



3.2 Main mathematical problems in TAT 

As in all tomographic methods, the following questions play the central 
role: 

Uniqueness of reconstruction: Is the collected data g sufficient 
for the unique reconstruction of the tomogram (function f{x))? 
When the author first looked at this problem in terms of the 
system ([l]), he was confused for a second. Indeed, the measured 
data g seems to be the boundary value of the solution of the 
wave equation in a cylinder, and the function / to recover is the 
initial condition. This clearly is impossible, since g essentially 
does not carry any information about /, besides the standard 
junction conditions where these functions meet. However, the 
thing is that ([1]) is not a boundary value problem in a cylinder, 
but rather a problem in the whole space, whose solution we ob- 
serve on the surface S only. What is even more important, if the 
sound speed is non-trapping (e.g., constant) and 5" is a closed 
surface (i.e., the boundary of a bounded domain), then the local 
energy decay theorems [35 |[5n[TT9lll37pi38j show that the solu- 
tion decays in the interior of S. This turns out to be the main 
feature that brings about uniqueness and inversion [2ll68p69pi24] . 

The question of sufficiency of the data collected on a non-closed 
hyper-surface S is much more complicated and can be considered 
largely resolved only in 2D, when the sound speed is constant 
and the function f{x) is compactly supported (see [6]i . as well 
as surveys [31I771I80] for the results and further discussion). This 
is a fascinating mathematical problem, which involves a variety 
of techniques from PDEs, integral, differential, and algebraic ge- 
ometry, microlocal analysis, commutative algebra, etc. A more 
general problem (with spheres replaced by the level surfaces of 
a polynomial) is discussed in L. Ehrenpreis' book [52]. 

Reconstruction procedures: Having a uniqueness of reconstruc- 
tion theorem is far from being sufficient for tomography, since 
one needs to be able to actually reconstruct the tomogram. So, 
the natural question, after proving uniqueness, is: how can one 
actually invert the mapping / i— )• 5? We will briefly describe the 
known algorithms in the next sub-section. 
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Some microlocal arguments in [B], although correct, were incomplete. The missing 
arguments are provided in [123j. 



Stability: Having uniqueness of reconstruction theorem, or even a 
reconstruction algorithm is no guarantee that the reconstruc- 
tion is practicably feasible (or at least that it can provide sharp 
pictures). The word "stability" in this context means the ef- 
fect that small errors in the measured data have on the recon- 
structed image. In other words, stability means well-posedness 
(in the Hadamard's sense j46j ) of the inverse problem. Unfortu- 
nately, essentially all tomographic problems are ill-posed to some 
degree. Some of them, such as Radon transform inversion, are 
very mildly unstable, which allows for wonderful CT-scan im- 
ages. Some others, like EIT and OT, are known to be severely 
ill-posed, and thus reconstruction of sharp images is practically 
impossible. 

The stability of the TAT/PAT reconstruction with full data (i.e., 
for the observation surface surrounding the object completely) 
is known to be very good, the same as for the Radon transform 
inversion, which leads to excellent reconstructions. This applies 
both to the case of an acoustically homogeneous object (when 
c{x) = const), as well as to the case of a variable non-trapping 
sound speed c{x) [2l[69l[l06l[l24l[l26l[l27] . The proofs are based 
upon inversion procedures and/or a microlocal analysis of the 
problem. 

Range: As it is common in integral geometry and tomography [GOpGl^ 



l63l[Ml[T03lll04j . the range of the forward operator f ^ g has in- 
finite co-dimension in the natural scales of function spaces. The 
description of this range is an important part of integral geo- 
metric and tomographic studies [ibid] . In the case of a spherical 
observation surface S and constant sound speed, the complete 
range descriptions are known [Il[4l[5l[9l[57j . Some of these range 
conditions are known to be necessary for more general obser- 
vation surfaces and sound speeds, but in these cases complete 
range description is not known (and might be impossible). 

Incomplete data: In the TAT/PAT case, one usually mentions an 
incomplete data situation, when the observation surface S 
does not surround the object completely [3l[77 1 [80 lll46U147j . This 
is a common, albeit somewhat misleading description, since de- 
pending on the location of the object, the "incomplete" data 
might be sufficient for unique, and sometimes even stable re- 
construction. It is quite common that even a small observation 
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surface S can provide enough data for proving uniqueness of 
reconstruction. For example, in the case of a constant sound 
speed, it is known that any set S that is not a part of an alge- 
braic hyper-surface, leads to the injectivity of the spherical mean 
operator Rg, and thus to unique TAT reconstruction (see [6] and 
references therein). Note that such S can be geometrically very 
small. It is clear that for all practical purposes, in spite of an 
uniqueness theorem, something should go wrong with the actual 
reconstruction in this case. And indeed, microlocal arguments 
similar to those of X-ray CT |118j show that most of the singu- 
larities (i.e., wave front set directions) of f{x) will be "invisible" 
("not audible") p[6m[77l[8m [Tn6l[TTT|[TT2l[TMll46[ll47j . This im- 
plies, in particular, the high instability of the reconstruction (all 
"invisible" interfaces will be blurred, no matter how sophisti- 
cated the algorithms are) [69 |[T06l[TM[Ti6l[Ti7] . This is also weh 
known in the more standard X-ray CT and SPECT |8 Hlll8j . 
The wonderful feature of the microlocal analysis is that it not 
only explains, but also allows one to predict these blurring effects 
[ibid]. 

On the other hand, if one is in the situation where there is a 
uniqueness result and all singularities of the object are "visible", 
one indeed can reconstruct the object stably [I3l[l4l[86]. Con- 
ditions of unique reconstruction and "visibility" have been also 
figured out for the case of variable sound speeds and can be ex- 
pressed simply in terms of geometric optics rays [3 ll69p77p80UlQ61 

While the questions above are common for all tomographic tech- 
niques, there is one issue that is specific for TAT (a similar, still not 
completely resolved, problem arises also in SPECT, see [76] and ref- 
erences therein, as well as the recent papers |20ll34j ). This is the 

Speed recovery problem: The problem ([1]) involves the unknown 
function /(x), the tomogram, as well as the sound speed c(x), 
which in all reconstruction methods is assumed to be known. 
It has been observed [691170] that errors in the values c(x) in- 
troduce significant artifacts into reconstruction. In other words, 
one needs to know the speed c{x) well, which is normally not the 
case. One of the options suggested to alleviate this difficulty, is 
to run an ultrasound transmission scan beforehand, which would 
provide the speed map [7^. However, there is numerical evidence 
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[148yi50j that it might be possible to determine both the speed 
c{x) and then the tomogram f{x) from the TAT data. Prov- 
ing uniqueness of reconstruction of c{x) happens to be a difficult 
problem. Some very limited partial results have been obtained 
recently (and not published, except [BTIES]) by various mathe- 
maticians: M. Agranovsky, D. Finch, K. Hickmann, Y. Hristova, 
P. Kuchment, L. Nguyen, P. Stefanov, and G. Uhlmann. How- 
ever, the problem (which, according to D. Finch's observation, is 
closely related to also not completely resolved well known trans- 
mission eigenvalue problem [36lB3lH5ll 7411 1 1 II 1 29 j ) is essentially 
open. Somewhat more understanding has been achieved concern- 
ing possible instability of the speed reconstruction (L. Nguyen, 
P. Stefanov, G. Uhlmann, unpublished). 

3.3 Reconstruction methods in TAT/PAT 

What techniques are available for actual TAT reconstructions? There 
are several groups of approaches that have been used successfully: 

Closed form inversion formulas: Closed-form inversion formulas 
in tomography bring about better theoretical understanding and 
lead to efficient reconstruction algorithms. The best known ex- 
ample is the so-called filtered backprojection (FBP) algorithm 
in X-ray tomography, which is derived from one of the popu- 
lar inversion formulas for the classical Radon transform (see, for 
example [63ti65l [T03l[TM] 1 . 

The first such formulas for TAT (the very existence of which had 
not been clear) were obtained in odd dimensions in [55], for the 
observation surface S being a sphere and the sound speed being 
constant. Let the function /(x) be supported within a ball of 
radius R and the detectors be located on the boundary S = dB 
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of this ball. Some of the formulas obtained in 1551 are: 



^(^) = "8^^^ y \y-x\ ""^^y^^ (') 



9B 



dB 



dA{y), (5) 

■r=|j/— x 



/(^O^-^/fi^f..*"'^'-' 



Stt'^R J \r dr \ dr 
dB 



dA{y), (6) 

= |j/-x| 



where dA{y) is the surface measure on 9-B and g represents the 
values of the spherical integrals ([3|) . 

Here differentiation with respect to r in ([5]) and ([6]) and the 
Laplace operator in (j3|) represent the filtration step, while the 
(weighted) integrals correspond to the backprojection: integra- 
tion over the set of spheres passing through the point x and 
centered on 5. 

A different family of inversion formulas valid in any dimension 
was found in [84] . Still another set of closed- form inversion for- 
mulas applicable in even dimensions was found in [54] • Finally, 
a unified family of inversion formulas was derived in |105j . 
Having closed form inversion formulas has the advantage that 
they usually lead to fast and precise inversion algorithms. How- 
ever, there are several disadvantages of these formulas in TAT. 
First of all, the FBP formulas are now available only for the 
observation surface being a sphere (see discussion above), with 
the only exception of [87], where such formulas are derived for 
a cube and some other crystallographic domains. Secondly, it is 
known (e.g., [77]) that if some part of the source function f{x) is 
supported outside the observation surface 5, then its reconstruc- 
tion inside S using FBP formulas might be incorrect. Finally, 
there are no FBP formulas known for the case of a variable sound 
speed. 

Eigenfunction expansions: - This approach, which theoretically 
works for arbitrary closed surfaces, was proposed in [85] (and 
extended in [2] to the case of variable sound speeds) . It is based 
on expansion into eigenfunctions of the Laplacian operator in the 
interior B oi S with zero Dirichlet conditions on 5. It is thus 
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nicely implementable whenever the spectrum and eigenfunctions 
of the Dirichlet Laplacian are known exphcitly, e.g., for spheres, 
half-spheres, cylinders, cubes and parallelepipeds, as well as the 
surfaces of some crystallographic domains. 

The function /(x) can be reconstructed inside B from the data 
g in ([T|), as the following L^(i?)-convergent series: 

f{x) = Y,fkMx)^ (7) 

k 

where tpk{x) are properly normalized eigenfunctions of the oper- 
ator (— c^(a;)A) in B with zero Dirichlet conditions and A^ are 
the corresponding eigenvalues. The Fourier coefficients fk can 
be recovered from the data 5 in ([3]) using one of the following 
formulas: 

00 
fk = XfgkiO) - Xf f sin {\kt)g'^{t)dt, 



00 

fk = XfgkiO) + Xf I cos {Xkt)g'k{t)dt, or 



00 00 



fk = -A^i /sin(Afct)<7fc(t)dt = -X-^ J J sm{\kt)g{x,t)^{x)dxdt, 
OS 

(8) 
where 

9k{t) = / g{x,t)—^{x)dx. 

s 

This method becomes computationally efficient when the eigen- 
values and eigenfunctions are known explicitly and a fast sum- 
mation formula for the series ([7]) is available, for instance when 
the acquisition surface 5 is a surface of a cube, and thus the 
eigenfunctions are products of sine functions. The resulting 3D 
reconstruction algorithm is extremely fast and precise (see |85j). 
This method applies in any dimension and is stable. It also 
does not have the deficiencies of the FBP formulas that we have 
mentioned above. Namely, presence of a part of the function 
f{x) outside f{x) does not hurt the reconstruction inside. The 
method, at least theoretically, works for arbitrary closed observa- 
tion surface s and variable speed c{x). However, its practicality 
in these circumstances is still questionable. 
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Time reversal: - We now describe an inversion technique that has 
the same advantages as the eigenfunction expansion method above, 
and in addition is easy to implement for any shape of the obser- 
vation surface and acoustically inhomogeneous media. One can 
come up easily with this method, if one notices the underlying 
assumption of TAT, which is often hidden, and which we have 
discussed explicitly before: local energy decay of the solution of 
dl]). Then one can naturally come up with the idea of running 
the wave equation in ([T|) back in time, starting at the infinite 
time with zero initial condition (which reflects the local energy 
decay) and using the measured data g as the boundary value. 
Eventually, at time t = 0, one arrives to the tomogram f{x). 
Certainly, the solution never vanishes exactly at any finite time 
(unless the sound speed is constant and the dimension is odd, 
where the Huygens' principle kicks in). Thus, one has to start 
from some sufficiently large time t = T with zero conditions and 
go back in time to arrive to an approximation of f{x), which one 
expects to get better when T — )■ oo. This approach has been im- 
plemented by various researchers, its feasibility was shown, and 
error estimates were provided (see, e.g., [62 t [68 l[69tll45| ). This 
works |62^l69j even in 2D (where decay is the slowest) and in 
inhomogeneous media. However, when trapping occurs, some 
parts become "invisible" and blur away |69j . 
A more sophisticated (than just zero) cut-off at time t = T is 
used in the version of time reversal suggested in |1171ll24tll25j . 
This leads to an equation with a contraction operator, which 
allows the use of the Neumann series (fixed point iterations) to 
obtain high quality images. 

It is the author's belief that the time reversal is the most versatile 
and easy to implement method of TAT reconstructions. 

Algebraic reconstruction techniques: ART is the well used workhorse 
in approaching inverse problems (especially, when analysis of the 
problem is too complicated). To put it simply, in ART one dis- 
cretizes the problem and uses one's favorite method (usually an 
iterative one) for solving the resulting linear algebraic system. 
Such techniques have been used in CT for quite a while [103j . 
ART algorithms frequently produce very good images. However, 
they are notoriously slow. In TAT, they have been used success- 
fully for reconstructions with partial data |13^lll3j and sound 
speed recovery [148yi50j . 
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Parametrix approaches: Some of the earher non-iterative recon- 
struction techniques of approximate nature |75 1I1151I116| were 
based (exphcitly, or implicitly) upon microlocal analysis. For 
example, in [75], by approximating the integration spheres by 
their tangent planes at the point of reconstruction and applying 
one of the inversion formulas for the classical Radon transform, 
one reconstructs a decent approximation to the image. 

Such techniques are related to the general scheme proposed in 
|29j for the inversion of the "generalized" Radon transform that 
integrates over curved manifolds. One constructs a parametrix 
(usually an integral Fourier operator (FIO)) for the forward op- 
erator F : f ^^ g, i.e. such operator P that the operators PF — I 
and FP — I, while not equal to zero, as in the case of a true in- 
version, are "smoothing." Thus, applying a parametrix P to the 
data g, one recovers the image / up to addition of a smooth func- 
tion. This also often reduces the problem to a Fredholm integral 
equation of the second kind, which is well amenable to numeri- 
cal solution. In other words, the parametrix method provides an 
efficient pre-conditioner for an iterative solver; the convergence 
of such iterations can be much faster than that of algebraic iter- 
ative methods. On the other hand, parametrix reconstructions 
can be often accepted as approximate images. 

3.4 Examples of TAT/PAT reconstructions 

Now, the mathematics looks nice, but does the TAT/PAT procedure 
really work? The reader can find below (Fig. [2]) the wonderful PAT 
image (curtesy of Wikipedia) of a blood vessel structure inside a hu- 
man hand. 

The next Fig. [3] (from |146j ) shows several TAT reconstructions of 
a muscle+fat phantom shown on top left. 

3.5 Quantitative PAT 

The usual TAT/PAT procedure essentially recovers the initial pressure 
f{x), which is proportional to the energy deposition function H{x). 
How does H{x) relate to the actual electric or optical parameters of the 
medium? It is less of a problem for TAT, where radio frequencies are 
used and one can hope to obtain more or less homogeneous irradiation 
of all tissues. The situation is different in PAT. In the diffusion regime. 
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Figure 2: A PAT reconstruction of the blood vessel system of a human hand. 
The picture is courtesy of Wikipedia. 
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Figure 3: a. The phantom, b. Reconstruction from partial data with "invis- 
ible" details blurred, c. Reconstruction from partial data with all features 
visible, d. Reconstruction from full data |146j . 

one has 

H{x) = r{x)cr{x)u{x), 

where F is the so called Griineisen coefficient, cr(x) is the EM en- 
ergy absorption coefficient, and u{x) is the radiation intensity. The 
following equation is satisfied by u{x): 

-V • D{x)S/u{x) + a{x)u{x) = 0, 
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where D is the diffusion coefficient. The question arises whether, after 
doing the TAT/PAT reconstruction and recovering H{x), one can go 
further and recover the actual optical parameters (-D, a, T) from H{x)l 
This is the goal of the so called Quantitative PAT (QPAT), which has 
started developing only recently [2T | [22l[2il[26llT7A9| [T20l[Tl8] and is 
in a very active stage now. In other words, QPAT takes of where PAT 
lands. 

4 Acousto-Electric Tomography 

Electrical impedance tomography (EIT) strives to recover the interior 
distribution of electric conductivity by measurements on the bound- 
ary [lllETlEIlETlillliaEallMlETlllilliaiiMlM Namely, one 
creates various boundary voltages and measures the resulting bound- 
ary currents (or vice versa). From these measurements one tries to 
recover the internal conductivity. The mathematical incarnation of 
EIT is the inverse conductivity problem, which was apparently 
suggested first by E. Calderon and has by now a glorious 30 years his- 
tory. Efforts of many leading mathematicians were directed towards 
proving that the measured data is sufficient for the unique recovery 
of internal conductivity. This happens to be much more mathemat- 
ically difficult topic than those arising in traditional tomography or 
in TAT. To large extend, this is due to the significant non-linearity 
of the corresponding mapping and instability of its inversion. This 
is still an active area, for instance since the optimal result in 2,D is 
still not known. However, the general understanding is that one does 
have sufficient information for the reconstruction of the conductiv- 
ity. We provide just a sample of references to this bursting with life 
topic [ISIEIIETIIIIIEHIIISOIIISSIIIMIIIM], where the reader can find 
plenty of information and further references. 

In this text we will be interested in the issues related to the actual 
reconstruction in EIT, which has also attracted enormous attention of 
scientists, including mathematicians. Due to the already mentioned 
instabilitjo, the pictures come out very much blurred and of low res- 
olution. This is a rigorously proven fact of life, not the deficiency of 
the mathematics used. It is the author's opinion that by now experts 
have achieved as good EIT reconstructions as humanly possible ... un- 



^Here instability = ill-posedncss ~ exponential decay of singular values of the direct 
mapping. 
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Figure 4: The AET procedure: electrical boundary measurements are done 
concurrently with scanning the object with ultrasound. The picture is cour- 
tesy of L. Kunyansky. 

less one changes the physical set-up of the measurements and/or uses 
some a 'priori information. Changing the experimental set-up is what 
is suggested in Acousto-Electric Tomography (AET). 

The main obstacle of EIT (and similarly of OT (optical tomogra- 
phy)) is that the signals measured at the boundary loose, exponentially 
fast with the depth, the information on where they came from. This is 
a simple minded explanation of why reconstructions are blurred. So, 
if one could somehow get some type of "interior" information about 
where the signals came from, one could hope for stabler image recon- 
structions. AET is one of the implementations of this idea. It has 
been known for some time that ultrasound irradiation of soft tissues 
modifies the tissues' electric and optical properties (electro-acoustic 
effect [88l[89]). It was thus an easy step to decide to send an ultra- 
sound beam that focuses on some internal location x and thus modifies 
(by a multiplicative factor close to 1) the electric conductivity (y{x) at 
this location. This would lead to a perturbation of the boundary EIT 
measurements, and what is crucial, the practitioner will know where 
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the perturbation came from - from the pomt x. Then one could scan 
the focused beam throughout the whole object and get hopefully suf- 
ficient information for a stable reconstruction. This idea of AET was 
suggested and tried by a direct measurement in |149j . It was shown 
there that a detectable (albeit rather small) signal does exist. How- 
ever, no reconstruction was done at that time. In a few years, the 
topic started developing fast [TOKIIlEOllMlESlESlEDKIlH], sometimes 
with the researchers being unaware of the original work |149| . Let us 
describe briefly the current state of affairs in AET (although by the 
time of publication, the situation will definitely change). 

The following observation was made experimentally and justified 
theoretically p8 | [89 l [149] : the acousto-electric effect, although de- 
tectable, is so small, that one can safely linearize the problem. 

Another smallness assumption was used in most works: the abil- 
ity of sharp focusing at a given location (i.e., creating a delta-type 
ultrasound pulse). Such perfect focusing is clearly impossible (see the 
discussion in the book [66] devoted to this issue). Still, let us assume 
for the time being that sufficiently good focusing is possible (and re- 
turn to this discussion in Section [5]) . This allows one to use the well 
studied "small volume inclusion" asymptotics, as in p^ [TT ] [38 l [39]. 
where such asymptotics play the crucial role. On the other hand, 
in [THllTS], only smallness of the acousto-electric effect is needed, and 
no perfect focusing is required (see Section [5]) . 

In all these works the authors studied what kind of interior quanti- 
ties can be stably recovered from the measurements, if perfect focusing 
(in particular, small volume asymptotics) is possible. For instance, if 
ui{x),U2{x) are the (unknown) potentials created by some boundary 
current set-up, then one can recover the values 

a{x)'Vui{x) ■ Vu2{x) 

for any interior point x, where u ■ v denotes the inner product of two 
vectors. Some other local functionals of the form F{a{x),Ui{x), Vui{x)) 
could be recovered (this is also the case in the previously mentioned 
MREIT and GDI pUTnnlfT22lfn3] . which we do not discuss in this 
paper). 

It was shown then that such values, if recovered from measure- 
ments, lead to locally unique and stable reconstruction of the conduc- 
tivity a [30 ] l38 U 78 U 79 yi28j . Essentially, one can prove that the Frechet 
derivative of the mapping 

a I—)- values of F{a, Ui, Vui) 
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(in appropriate function spaces) is an injective semi-Fredholm opera- 
tor. 

A variety of inversion procedures has been suggested and mostly 
tested on numerical phantoms: those involving numerical optimization 
[TT1[38], those reduced to solving well posed hyperbolic problems [16], 
or the ones that lead to solving transport equation or Poisson type 
elliptic equations [781179]. 

In most cases one could achieve wonderful quality reconstructions, 
e.g. see Fig. \5\ where the method of [79j is used, which reduces to 
solving a Poisson equation for determining the conductivity. 
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Figure 5: An example of AET reconstruction. The picture is courtesy of L. 
Kunyansky. 

Looking at Fig. [5l one might (and should) feel cheated. Indeed, 
how can one get such a good reconstruction with the data contami- 
nated by a 50% noise? The answer will be given in the next section. 

By now, the mathematics of AET, albeit just a few years old, is 
already rather successful. The experimental implementation of the 
AET lags behind, due to the difficulty of acquiring good signal-to- 
noise ratios. 

A different combination of ultrasound and EIT is suggested in [59] . 
Here one again creates currents through the interior of the body of in- 
terest. These currents lead to a small inhomogeneous heating of the 
tissues, and thus to thermoelastic expansion. Then the TAT proce- 
dure, using the microphones surrounding the body, reconstructs a lo- 
cal functional F{a{x),u{x),'Vu{x)), after which one of the previously 
mentioned procedures of reconstruction can be applied. 
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5 Synthetic focusing in hybrid tech- 
niques 

As it has been mentioned, the unfeasible [66] perfect focusing is as- 
sumed in most mathematical work on AET, in particular, when using 
the small volume asymptotics. Can this be avoided? The answer, as 
it was explained in [75] and then confirmed in [73], is a "yes." Indeed, 
the delta functions that are idealized focused beams, form a function 
"basis." Suppose that we can produce another, practically feasible 
set of ultrasound waves, which would also form such a basis. Then, 
using the smallness of the acousto-electric effect, and thus linearity, 
one could mathematically process the data obtained from that ba- 
sis of waves and "synthetically focus" them by changing the basis to 
delta- functions. 

There are several examples of possible bases from [78j : 

• Using large planar broad band transducers, one could generate 
a set of monochromatic planar waves with arbitrary wave vec- 
tors, and then the synthetic focusing would be just applying the 
inverse Fourier transform. This option was adopted in |16j . Its 
practical feasibility is not yet clear. 

• Using point-like omni-directional transducers, one could generate 
thin spherical shell waves. Then, lo and behold, the synthetic 
focusing will boil down to inversion of a restricted spherical mean 
transform, and thus any of the standard TAT inversions would 
do it. This is the option of [79]. The problem with this is that 
it is much easier to create a short A^-shaped (Fig. [6|) spherical 
wave rather than (5-shaped such wave. 




Figure 6: An A^-shaped pulse. 

• One can create narrow "pencil beams", as it is done in [142j . 
and then synthetic focusing would coincide with the inversion of 
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the standard X-ray transform. This option has to struggle with 
impossibihty of creating a truly homogeneous pencil ultrasound 
beam of sufficient length. [66j 

The sharp and amazingly stable with respect to the noise recon- 
struction shown in Fig. [5l was done using the second option for 
the basis: thin spherical waves, with consecutive TAT inversion for 
the synthetic focusing. After that an elliptic (Poisson type) equation 
was solved to recover the conductivity. Now, what about having A^- 
shaped rather than (5-shaped pulses? This "difficulty" turns out to be 
a blessing. Indeed, the TAT reconstruction includes a filtering portion, 
which increases the noise and decreases stability. However, with the 
A^-shaped pulse, this filtration step is not needed, since it is already 
performed by the transducer. As the result, the synthetic focusing by 
the TAT inversion happens to be a smoothing operator and thus 
kills a lot of noise. If we could indeed produce and use (5-shaped pulses, 
the reconstruction would work, but would be very unlikely to survive 
a 50% noise. 

6 Ultrasound Modulated Optical To- 
mography (UMOT) 

The idea of scanning an object with a focused ultrasound that we 
applied in AET, can be tried with the optical tomography (OT) as 
well. The goal is the same: to improve drastically the resolution of 
OT (which is dismal on a centimeter depth and deeper). Since OT, 
like EIT, is a cheap, safe, and high contrast modality, achieving this 
goal, and thus adding high resolution, would make it an invaluable 
diagnostic tool. 

In comparison with the AET, the situation with UMOT is re- 
versed: there is an extensive body of experimental research (see |142j 
and references therein), but the first glimpses on the mathematics of 
UMOT are just appearing [8 | l25 yiOHI102j . and even the mathematical 
model is not settled down. 

The set-up of OT is as follows: One sends a beam of (coherent 
or incoherent) laser light through the body of interest and observes 
the intensity and speckle patterns of the outgoing light. The features 
of interest are the internal distribution of the absorption and scatter- 
ing coefficients. The contrast in optical properties of cancerous and 
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healthy locations is often huge. However, diffused photons, when they 
reach a detector after multiple scattering essentially loose any infor- 
mation about the locations they went through. This leads to dismal 
resolution at centimeters' depth (although good pictures can be ob- 
tained at skin depth). The idea of UMOT is to, concurrently with 
OT measurements, scan the body with a focused ultrasound and thus 
to acquire some interior (i.e., location-dependent) information, which 
hopefully would stabilize the problem. It has been argued in physics 
and biomedical engineering literature (e.g., |72 t l90 1[142j ) that when 
using coherent light and measuring the ultrasound frequency Fourier 
component of the outgoing speckle pattern, one can recover the values 
of the following functional at an interior location x: 

G{x,d)A'^{x)I{x). 

Here A is the applied ultrasound power (assumed to be known), I is 
the light intensity, d is the detector position on the boundary, and 
G{x, d) is the "probability of a photon emitted at the location x to 
reach the detector at the location d". In other words, G{x,d) is a 
Green's function of the diffusion equation 

-VD{x)VI{x) + ^ia{x)I{x) = 

inside the domain of the interest. The difficulty (at least for the au- 
thor) is determining what the "probability of reaching the detector" 
means (e.g., does this mean the first time of reaching the point dl). 
Thus, it is not clear what boundary conditions the Green's function 
should correspond to. 

It was assumed in [8] that the correct boundary conditions are 
those that correspond to the optical impedance at the boundary of the 
object. Under this condition and with the perfect focusing assump- 
tion, a reconstruction algorithm was applied that showed sufficiently 
sharp internal reconstructions of the absorption coefficient /x^ (al- 
though the quality was lower than in the AET case) . It was also shown 
in [8] (see also the acknowledgments there) that (formally computed) 
Frechet derivative of the forward mapping is a semi-Fredholm opera- 
tor in natural function spaces. However, injectivity of this derivative 
was not shown. Thus there are so far no local injectivity results. 

Some controversy surrounds the usage of coherent light. It is 
claimed in engineering literature |142| that the signal from ultrasound 
modulation in the case of incoherent light could not be detected so 
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far. However, there already are some mathematical studies of UMOT 
using incoherent light |25j . 

Synthetic focusing in UMOT is possible. However, while the spher- 
ical waves should still do the job, the use of planar waves and the con- 
sequent inversion of the Fourier transform seems to be not an option 
here, due to the presence of the square of the acoustic power A{x) in 
the measured functional. 

7 Why the improvement? Inverse prob- 
lems with interior information 

The examples of AET and UMOT show how acquiring interior (i.e., 
attached to the internal locations) information stabilizes the utterly 
unstable inverse problems. This is an example of a folklore meta- 
statement: "appropriate" internal information stabilizes the severely 
unstable problems like diffused OT or EIT. 

This issue was studied in detail in [16] , where several internal infor- 
mation functionals arising in applications (including those described 
above) were studied. Different functionals required sometimes differ- 
ent techniques. There is a feeling though that there might be an an- 
swer to a general question: What kind of a function F{D{x), a{x),u{x), Vu{x)), 
if known, stabilizes the inverse boundary problems for 

-V ■ D{x)Vu + au = 07 

A rather general answer was just given in [128] . Under some rea- 
sonable conditions on the functional, which cover all cases we have 
considered here, it is shown that the Frechet derivative of the for- 
ward mapping is a semi-Fredholm operator. This explains why the 
observed improvement in stability occurs. Together with local unique- 
ness (which probably need to be proven individually in each case) this 
also gives local uniqueness and stability. 
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